Este sitio web aborda el análisis de sólidos suspendidos totales (SST) en cuerpos de agua mediante el uso de imágenes satelitales y técnicas de aprendizaje automático. Presenta una introducción teórica sobre la importancia de los SST como indicador ambiental y luego desarrolla una parte práctica con Python, donde se integran datos espectrales y mediciones reales para aplicar modelos de regresión. Se evalúa el desempeño de estos modelos y se exploran relaciones entre bandas espectrales y SST, proponiendo mejoras metodológicas para estudios futuros.
Palabras clave
GISTAQ, UTN, FRRe, Quarto, Sólidos Suspendidos
1 Sólidos suspendidos totales, por Vera Geneyer (https://github.com/VeraGeneyer)
Los sólidos suspendidos totales (TSM) son la cantidad de materia en suspensión en el agua, que incluye plancton, minerales, arena, y microorganismos. Se determinan como el residuo no filtrable de una muestra de agua. Niveles altos (TSM) pueden reducir la transparencia del agua, limitar la luz y y transportar sustancias tóxicas, afectando la vida acuática y la calidad del agua. Este parámetro, medido mediante sensores remotos, nos da información sobre el estado físico del cuerpo de agua y están relacionados con factores como la humedad, temperatura y entre otros, que es vital para detectar riesgos al ecosistema y cumplir con las normas ambientales.
1.1 Métodos tradicionales
Tabla 1: Características principales de algoritmos tradicionales para la estimación de sólidos suspendidos.
De acuerdo a un estudio que analizó 48 cuerpos de agua, la estimación de TSM se hizo en su mayoría por modelos lineales, siendo la banda B8A la más frecuente [3].
1.2 Métodos de aprendizaje automático
El aprendizaje automático (ML) es una rama de la inteligencia artificial cuyo objetivo es desarrollar algoritmos capaces de resolver problemas mediante el análisis de datos y la creación de funciones que describen el comportamiento de fenómenos monitoreados [4]. Los modelos de aprendizaje automático más utilizados y mencionados por los investigadores para predecir la concentración de SST son:
Bosque Aleatorio (RF) y Refuerzo Adaptativo (AdB), modelos que se destacan por su robustez ante datos complejos y ruidosos. Estos algoritmos construyen múltiples árboles de decisión que analizan las relaciones entre características como el uso del suelo o el volumen de escorrentía y los niveles de SST [5].
Redes Neuronales Artificiales (ANN), copian las redes neuronales biológicas y aprenden patrones complejos en grandes volúmenes de datos, como los niveles de SST en distintas condiciones ambientales [5],
k-Nearest Neighbors (kNN), en sus variantes de ponderación uniforme y variable, que estima el SST en función de la cercanía en características de nuevos puntos de muestreo con datos históricos [5].
El aprendizaje automático es esencial para mejorar la precisión y rapidez en el análisis de la calidad del agua, proporcionando un monitoreo más eficiente y menos costoso en comparación con los métodos tradicionales, especialmente en áreas de difícil acceso o con datos limitados.
Tabla 2: Características principales de algoritmos de aprendizaje automático para la estimación de sólidos suspendidos.
Modelo de machine learning
Software
Agua
Datos
Métricas
Referencias
Bagging y Random Forest
Programa R
Bahía
Muestreo
Prueba de normalidad multivalente Mardia-tests y Royston
Regresión lineal, LASSO, regresión de vectores de soporte (SVR), K vecinos más cercanos (KNN), bosque aleatorio (RF) y redes neuronales artificiales (ANN).
En esta sección describo el proceso completo para crear un modelo que estima los sólidos suspendidos totales (SST) usando datos satelitales y de campo, desde la preparación de datos hasta la validación del algoritmo.
2.1 Procesamiento inicial de datos con Sen2Cor
Se detalla el procedimiento técnico que implementé para procesar información ambiental georreferenciada, con el objetivo de analizar el comportamiento del parámetro sólidos suspendidos (sol_sus) en una región específica del estudio (píxel 3x3).
Para esto, utilicé el lenguaje Python y la biblioteca pandas, que resulta particularmente eficiente para el manejo y transformación de datos en estructuras tabulares. Este paso fue fundamental para preparar los datos y asegurar su calidad antes del análisis.
2.1.1 Carga de datos
Primero importo la biblioteca pandas, una herramienta en Python que se utiliza para manejar datos en formato tabular (como hojas de cálculo o CSVs). Se le da el alias pd por convención, para simplificar el código.
Luego cargo dos archivos CSV con la función pd.read_csv(), la cual convierte dichos archivos en objetos del tipo DataFrame, que representan tablas en memoria, que son estructuras de datos similares a tablas (parecida a una hoja de Excel). Los conjuntos de datos cargados fueron:
lab_df: contiene datos de laboratorio, incluyendo el parámetro de interés sol_sus.
Verifico la carga correcta mostrando las primeras filas con la función .head(). Es útil para ver rápidamente cómo es la estructura del archivo: qué columnas hay, qué tipo de datos, si se cargó bien.
import pandas as pd # pandas es la biblioteca para manejar datos tabulares# Cargar los archivos de datosgis_df = pd.read_csv('datos/base_de_datos_gis.csv')lab_df = pd.read_csv('datos/base_de_datos_lab.csv')# Ver las primeras filas para asegurarse de que se cargaron biengis_df.head(), lab_df.head()print("Primeras filas de gis_df:")display(gis_df.head())print("\nPrimeras filas de lab_df:")display(lab_df.head())
Primeras filas de gis_df:
fecha
punto
pixel
banda
reflect
longitud
latitud
0
2023-05-11
1
1x1
B01
0.161744
-58.868047
-27.464687
1
2023-05-11
1
1x1
B02
0.170400
-58.868047
-27.464687
2
2023-05-11
1
1x1
B03
0.208800
-58.868047
-27.464687
3
2023-05-11
1
1x1
B04
0.251000
-58.868047
-27.464687
4
2023-05-11
1
1x1
B05
0.257656
-58.868047
-27.464687
Primeras filas de lab_df:
fecha
longitud
latitud
param
valor
0
2023-05-11
-58.868047
-27.464687
ph
6.8
1
2023-05-11
-58.868047
-27.464687
cond
141.1
2
2023-05-11
-58.868047
-27.464687
sol_sus
198.0
3
2023-05-11
-58.868047
-27.464687
turb
185.5
4
2023-05-11
-58.866729
-27.466303
ph
6.8
📄 Nota técnica
pd.read_csv() carga los archivos en estructuras llamadas dataframes, que funcionan como tablas.
head() te muestra las primeras 5 filas para ver cómo están los datos.
display() permite mostrar las tablas con formato más visual (en HTML).
2.1.2 Filtrado del parámetro sólidos suspendidos (‘sol_sus’)
En el conjunto de datos del laboratorio lab_df, hay múltiples parámetros medidos (como pH, turbidez, etc.). En este caso, me interesa trabajar solamente con los datos de sólidos suspendidos, identificado como "sol_sus" en la columna param. Este filtrado selectivo lo realicé para limitar el análisis al fenómeno físico-químico de interés.
Filtré el DataFrame para quedarme solo con esas filas, y renombré la columna valor como sol_sus para que sea más claro en los siguientes pasos.
# Filtrar solo las filas donde el parámetro es 'sol_sus'sol_sus_df = lab_df[lab_df["param"] =="sol_sus"]# Renombrar la columna 'valor' a 'sol_sus' para que tenga sentido en el mergesol_sus_df = sol_sus_df.rename(columns={"valor": "sol_sus"})# Mostrar para confirmarsol_sus_df.head()print("Primeras filas de sol_sus_df:")display(sol_sus_df.head())
Primeras filas de sol_sus_df:
fecha
longitud
latitud
param
sol_sus
2
2023-05-11
-58.868047
-27.464687
sol_sus
198.0
6
2023-05-11
-58.866729
-27.466303
sol_sus
150.0
10
2023-05-11
-58.864889
-27.468056
sol_sus
101.0
14
2023-05-11
-58.863367
-27.469240
sol_sus
95.0
18
2023-05-11
-58.862454
-27.470561
sol_sus
69.0
📄 Nota técnica
lab_df[lab_df["param"] == "sol_sus"] filtra las filas cuyo valor en la columna "param" sea "sol_sus".
.rename(columns={"valor": "sol_sus"}) cambia el nombre de la columna "valor" a "sol_sus".
2.1.3 Reorganización de datos: pivotar bandas en columnas individuales
En este paso, convierto los valores únicos de la columna banda (como B01, B02, etc.) en nombres de columnas. Cada nueva columna contendrá los valores del parámetro reflect correspondientes a esa banda en particular. Esta operación se realiza antes de unir con los valores de sol_sus, ya que el valor de reflectancia depende de la banda, mientras que sol_sus es un dato independiente que se asignará luego por punto, fecha y ubicación.
# Pivotear la tabla para que cada banda sea una columnagis_pivot = gis_df.pivot_table( index=['fecha', 'punto', 'pixel', 'latitud', 'longitud'], columns='banda', values='reflect').reset_index()# Eliminar el nombre del índice de columnas generado por el pivotgis_pivot.columns.name =Noneprint("Primeras filas de gis_pivot:")display(gis_pivot.head())
Primeras filas de gis_pivot:
fecha
punto
pixel
latitud
longitud
B01
B02
B03
B04
B05
B06
B07
B08
B11
B12
B8A
0
2023-05-11
1
1x1
-27.464687
-58.868047
0.161744
0.170400
0.208800
0.251000
0.257656
0.202175
0.207500
0.194000
0.118625
0.114738
0.173069
1
2023-05-11
1
3x3
-27.464687
-58.868047
0.161290
0.171256
0.206644
0.250911
0.258058
0.202170
0.207490
0.193722
0.118892
0.115462
0.173137
2
2023-05-11
2
1x1
-27.466303
-58.866729
0.149925
0.170600
0.202100
0.243200
0.247169
0.186631
0.187512
0.175800
0.114100
0.111031
0.159263
3
2023-05-11
2
3x3
-27.466303
-58.866729
0.149906
0.170378
0.202278
0.242689
0.247115
0.186739
0.188027
0.175544
0.113992
0.111210
0.159373
4
2023-05-11
3
1x1
-27.468056
-58.864889
0.149725
0.166200
0.197100
0.223800
0.218175
0.159169
0.163344
0.150700
0.114263
0.113006
0.139662
📄 Nota técnica
pivot_table() reorganiza el DataFrame convirtiendo los valores de una columna (banda) en columnas individuales.
index=[...] define las columnas que se mantendrán como claves (se repetirán por fila).
columns='banda' indica qué columna queremos transformar en nombres de columnas.
values='reflect' especifica qué valor colocar en cada celda de la nueva tabla.
reset_index() convierte los índices jerárquicos en columnas normales para facilitar el análisis.
columns.name = None quita la etiqueta “banda” que se agregaría al encabezado por defecto.
2.1.4 Unión de datos geoespaciales con datos de laboratorio
Una vez que las bandas han sido transformadas en columnas, combino esta tabla con los valores de sólidos suspendidos (sol_sus) provenientes del laboratorio. La combinación se hace usando las columnas fecha, latitud y longitud, que permiten identificar los datos correspondientes a un mismo punto geográfico y temporal.
# Realizar el merge por ubicación y fechadf_merged = pd.merge( gis_pivot, sol_sus_df[['fecha', 'latitud', 'longitud', 'sol_sus']], on=['fecha', 'latitud', 'longitud'], how='inner')print("Primeras filas del DataFrame combinado:")display(df_merged.head())
Primeras filas del DataFrame combinado:
fecha
punto
pixel
latitud
longitud
B01
B02
B03
B04
B05
B06
B07
B08
B11
B12
B8A
sol_sus
0
2023-05-11
1
1x1
-27.464687
-58.868047
0.161744
0.170400
0.208800
0.251000
0.257656
0.202175
0.207500
0.194000
0.118625
0.114738
0.173069
198.0
1
2023-05-11
1
3x3
-27.464687
-58.868047
0.161290
0.171256
0.206644
0.250911
0.258058
0.202170
0.207490
0.193722
0.118892
0.115462
0.173137
198.0
2
2023-05-11
2
1x1
-27.466303
-58.866729
0.149925
0.170600
0.202100
0.243200
0.247169
0.186631
0.187512
0.175800
0.114100
0.111031
0.159263
150.0
3
2023-05-11
2
3x3
-27.466303
-58.866729
0.149906
0.170378
0.202278
0.242689
0.247115
0.186739
0.188027
0.175544
0.113992
0.111210
0.159373
150.0
4
2023-05-11
3
1x1
-27.468056
-58.864889
0.149725
0.166200
0.197100
0.223800
0.218175
0.159169
0.163344
0.150700
0.114263
0.113006
0.139662
101.0
📄 Nota técnica
pd.merge() permite combinar dos DataFrames en uno nuevo, uniendo filas que coincidan en las columnas especificadas.
on=["latitud", "longitud"] indica que la combinación debe hacerse usando esas columnas como claves.
how="inner" especifica el tipo de combinación:
"inner": solo conserva las filas donde hay coincidencia en ambos DataFrames.
Otras opciones:
"left": conserva todas las filas del primer DataFrame.
"right": conserva todas las filas del segundo.
"outer": conserva todo, incluso si no hay coincidencia.
2.1.5 Filtrado espacial por píxel de interés (3x3)
Luego de combinar los datos, aplico un filtrado adicional al DataFrame sobre la columna pixel para conservar únicamente las filas correspondientes al área geográfica designada como "3x3". Este paso reduce el dominio de análisis y permite concentrarse en una región de estudio concreta.
# Filtrar solo los datos del pixel 3x3df_pixel_3x3 = df_merged[df_merged["pixel"] =="3x3"]print("Primeras filas del pixel 3x3:")display(df_pixel_3x3.head())
Primeras filas del pixel 3x3:
fecha
punto
pixel
latitud
longitud
B01
B02
B03
B04
B05
B06
B07
B08
B11
B12
B8A
sol_sus
1
2023-05-11
1
3x3
-27.464687
-58.868047
0.161290
0.171256
0.206644
0.250911
0.258058
0.202170
0.207490
0.193722
0.118892
0.115462
0.173137
198.0
3
2023-05-11
2
3x3
-27.466303
-58.866729
0.149906
0.170378
0.202278
0.242689
0.247115
0.186739
0.188027
0.175544
0.113992
0.111210
0.159373
150.0
5
2023-05-11
3
3x3
-27.468056
-58.864889
0.149718
0.167256
0.197322
0.223511
0.218608
0.159292
0.163228
0.151544
0.114241
0.112931
0.139783
101.0
7
2023-05-11
4
3x3
-27.469240
-58.863367
0.148441
0.163289
0.188733
0.196167
0.183854
0.139388
0.140017
0.135400
0.115602
0.113097
0.128222
95.0
9
2023-05-11
5
3x3
-27.470561
-58.862454
0.145838
0.157733
0.178767
0.173911
0.162600
0.132140
0.133229
0.126656
0.117069
0.113992
0.124548
69.0
📄 Nota técnica
Usa filtrado booleano (DataFrame[condición]), que es la forma estándar en pandas para seleccionar subconjuntos de datos.
df_pixel_3x3 = df_combinado[df_combinado["pixel"] == "3x3"] selecciona ese subconjunto. Filtra las filas cuyo valor en la columna "pixel" es igual a "3x3".
2.1.6 Exportación de datos procesados
Finalmente, guardo el resultado como un nuevo archivo .csv dentro de la carpeta datos.
Por último, exporto el resultado a un nuevo archivo en formato .csv, mediante la función to_csv() de pandas, con el parámetro index=False para evitar que la columna de índice se incluya en el archivo de salida que pandas crea por defecto. Esto me permite utilizarlo después para visualización o análisis posterior.
# Guardar el archivo CSV dentro de la carpeta "datos"df_pixel_3x3.to_csv('datos/datos_sol_sus_pixel_3x3.csv', index=False)
📄 Nota técnica
to_csv() guarda los datos en formato CSV.
index=False evita que se guarde el índice numérico del DataFrame como una columna adicional en el CSV.
2.2 Aplicación del procesamiento a datos ACOLITE
A continuación, repetí el mismo procedimiento de procesamiento aplicado anteriormente, esta vez utilizando el archivo base_de_datos_gis_acolite.csv, generado con el procesador atmosférico ACOLITE.
ACOLITE es una herramienta diseñada específicamente para corregir efectos atmosféricos en ambientes acuáticos, a partir de imágenes satelitales. Gracias a que la estructura del archivo es similar al utilizado anteriormente (Sen2Cor), fue posible aplicar la misma lógica de filtrado, transformación y combinación para preparar los datos de reflectancia.
Código
import pandas as pd# 1. Cargar archivos CSV gis_acol_df = pd.read_csv('datos/base_de_datos_gis_acolite.csv')lab_df = pd.read_csv('datos/base_de_datos_lab.csv')# 2. Filtrar solo el parámetro 'sol_sus' sol_sus_df = lab_df[lab_df["param"] =="sol_sus"]sol_sus_df = sol_sus_df.rename(columns={"valor": "sol_sus"})# 3. Pivotear bandas en columnas gis_pivot = gis_acol_df.pivot_table( index=['fecha', 'punto', 'latitud', 'longitud'], columns='banda', values='reflect').reset_index()gis_pivot.columns.name =None# Quitar la etiqueta 'banda'# 4. Combinar reflectancias con sol_sus df_merged = pd.merge( gis_pivot, sol_sus_df[['fecha', 'latitud', 'longitud', 'sol_sus']], on=['fecha', 'latitud', 'longitud'], how='inner')# 5. Guardar resultado finaldf_merged.to_csv('datos/datos_sol_sus_acolite.csv', index=False)# 6. Mostrar muestra del resultadoprint("Primeras filas del DataFrame resultante:")display(df_merged.head())
Primeras filas del DataFrame resultante:
fecha
punto
latitud
longitud
B01
B02
B03
B04
B05
B06
B07
B08
B11
B12
B8A
sol_sus
0
2023-05-11
1
-27.464687
-58.868047
0.041367
0.046669
0.085630
0.128799
0.128704
0.091748
0.098365
0.083253
0.002521
0.000000
0.066612
198.0
1
2023-05-11
2
-27.466303
-58.866729
0.032235
0.047916
0.084191
0.123898
0.122481
0.079922
0.082403
0.069021
0.001669
0.000000
0.055924
150.0
2
2023-05-11
3
-27.468056
-58.864889
0.030956
0.044084
0.078444
0.104498
0.095737
0.052162
0.057475
0.044354
0.000209
0.000174
0.034767
101.0
3
2023-05-11
4
-27.469240
-58.863367
0.029646
0.040407
0.070384
0.078357
0.063937
0.032452
0.033994
0.028441
0.001279
0.000000
0.022597
95.0
4
2023-05-11
5
-27.470561
-58.862454
0.026737
0.034968
0.060701
0.056602
0.044067
0.024809
0.026684
0.019209
0.001879
0.000000
0.017930
69.0
2.3 Modelo base · Regresión lineal simple
En este análisis aplico un modelo de regresión lineal simple para estudiar la relación entre la reflectancia (banda B01) y la concentración de sólidos suspendidos (sol_sus), utilizando datos experimentales. La regresión lineal es una técnica fundamental del aprendizaje automático supervisado que permite predecir un valor continuo a partir de una o más variables independientes, y suele emplearse como modelo base frente a enfoques más complejos. A lo largo de esta sección, detallo paso a paso las acciones realizadas y explico los conceptos clave para comprender y replicar este análisis.
2.3.1 Importar librerías
En este paso, cargo las bibliotecas necesarias para procesar datos, ajustar modelos de regresión, evaluar su desempeño y visualizar los resultados.
import pandas as pdfrom sklearn.model_selection import train_test_splitfrom sklearn.linear_model import LinearRegressionfrom sklearn.metrics import mean_squared_error, r2_scoreimport matplotlib.pyplot as plt
📄 Nota técnica
pandas se utiliza para manejar datos en forma de tablas (DataFrames), especialmente útiles al trabajar con archivos .csv.
train_test_split permite dividir los datos en subconjuntos de entrenamiento y prueba, lo cual es esencial para evaluar el desempeño de un modelo sin sobreajustarlo.
LinearRegression representa un modelo lineal que se ajusta a los datos minimizando el error cuadrático entre las predicciones y los valores reales.
mean_squared_error y r2_score son métricas de evaluación: el primero mide el promedio de los errores al cuadrado, mientras que el segundo indica qué tan bien el modelo explica la variabilidad de los datos.
matplotlib.pyplot se utiliza para crear gráficos. Permite visualizar los datos y los resultados del modelo.
2.3.2 Cargar datos desde un CSV
Importo el archivo .csv con los datos experimentales. Se visualizan las primeras filas para verificar que los datos se han cargado correctamente.
# Cargar el CSVdatos = pd.read_csv('datos/datos_sol_sus_acolite.csv')# Mostrar las primeras filas para verificardatos.head()print("Primeras filas de datos:")display(datos.head())
Primeras filas de datos:
fecha
punto
latitud
longitud
B01
B02
B03
B04
B05
B06
B07
B08
B11
B12
B8A
sol_sus
0
2023-05-11
1
-27.464687
-58.868047
0.041367
0.046669
0.085630
0.128799
0.128704
0.091748
0.098365
0.083253
0.002521
0.000000
0.066612
198.0
1
2023-05-11
2
-27.466303
-58.866729
0.032235
0.047916
0.084191
0.123898
0.122481
0.079922
0.082403
0.069021
0.001669
0.000000
0.055924
150.0
2
2023-05-11
3
-27.468056
-58.864889
0.030956
0.044084
0.078444
0.104498
0.095737
0.052162
0.057475
0.044354
0.000209
0.000174
0.034767
101.0
3
2023-05-11
4
-27.469240
-58.863367
0.029646
0.040407
0.070384
0.078357
0.063937
0.032452
0.033994
0.028441
0.001279
0.000000
0.022597
95.0
4
2023-05-11
5
-27.470561
-58.862454
0.026737
0.034968
0.060701
0.056602
0.044067
0.024809
0.026684
0.019209
0.001879
0.000000
0.017930
69.0
📄 Nota técnica
pd.read_csv carga datos desde un archivo .csv y los convierte en un DataFrame de Pandas. Esta estructura tabular permite filtrar, seleccionar y transformar fácilmente los datos.
datos.head() permite ver las primeras 5 filas del DataFrame para tener una vista preliminar de los datos cargados.
2.3.3 Definir variables y particionar el conjunto
Selecciono las variables relevantes: B01 como variable independiente y sol_sus como variable dependiente. Luego divido el conjunto en dos subconjuntos: uno para entrenar el modelo y otro para probarlo, lo cual sirve para evaluar su capacidad de generalización.
# Selección de variablesX = datos[["B01"]] # variable independientey = datos["sol_sus"] # variable dependiente# División en entrenamiento y validaciónX_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=False)
📄 Nota técnica
Se selecciona una columna como variable independiente (X) y otra como variable dependiente (y). Es importante usar doble corchete al seleccionar una sola columna como X para mantener la estructura de tabla.
train_test_split divide el conjunto de datos en entrenamiento y prueba. Esto permite entrenar el modelo en un subconjunto y evaluar su capacidad de generalización con otro.
El parámetro test_size=0.2 indica que el 20% de los datos se usan para prueba. shuffle=False mantiene el orden original de los datos, útil cuando los datos están organizados temporalmente o espacialmente.
2.3.4 Entrenar modelo de regresión lineal
En este paso se entrena un modelo de regresión lineal usando los datos de entrenamiento. El modelo aprende la relación matemática entre la reflectancia y los sólidos suspendidos.
LinearRegression().fit() ajusta un modelo lineal a los datos de entrenamiento. Internamente calcula la pendiente e intercepto que minimizan la diferencia entre las predicciones y los valores reales.
2.3.5 Evaluar desempeño del modelo
Una vez entrenado el modelo, evaluo su desempeño usando métricas estadísticas. Estas nos permiten cuantificar qué tan bien el modelo predice los valores de sólidos suspendidos a partir de la reflectancia en los datos de prueba.
El error cuadrático medio es: 1846.79
El coeficiente de determinación (R²) es: -0.252
📄 Nota técnica
predict() genera predicciones del modelo usando los datos de prueba. Estas predicciones se comparan con los valores reales para evaluar el desempeño.
mean_squared_error calcula el promedio de los errores al cuadrado. Cuanto menor sea este valor, mejor se ajusta el modelo.
r2_score mide qué proporción de la variabilidad en los datos es explicada por el modelo. Un valor cercano a 1 indica una buena predicción.
2.3.6 Visualizar resultados
Finalmente, se visualiza gráficamente la relación entre reflectancia y sólidos suspendidos, tanto en el conjunto de entrenamiento como en el de prueba. Esto ayuda a interpretar de forma visual cómo se ajusta el modelo a los datos reales.
plt.subplots crea una figura con uno o más ejes para dibujar. Permite organizar varios gráficos en una misma figura.
plot() traza una línea continua. Se usa para mostrar la línea de regresión generada por el modelo.
scatter() traza puntos individuales. Se usa para mostrar los datos reales y compararlos con la línea del modelo.
set() configura etiquetas de ejes y títulos de los subgráficos.
legend() muestra una leyenda que identifica cada elemento del gráfico.
fig.suptitle() agrega un título general a la figura completa.
plt.show() es necesario para visualizar los gráficos al renderizar el documento.
2.4 Regresión lineal individual por banda
Con el fin de profundizar el análisis, evalué la relación entre los sólidos suspendidos (sol_sus) y cada banda espectral por separado. Para esto, entrené un modelo de regresión lineal simple usando los mismos datos experimentales para cada banda.
Con este enfoque comparo la capacidad predictiva individual de cada banda mediante las métricas R², R² ajustado y el error cuadrático medio (RMSE). Este análisis me permite identificar qué bandas tienen mejor relación lineal con los sólidos suspendidos y sentar la base para probar modelos más complejos.
Código
import pandas as pdimport numpy as npfrom sklearn.model_selection import train_test_splitfrom sklearn.linear_model import LinearRegressionfrom sklearn.metrics import mean_squared_errorimport matplotlib.pyplot as pltimport mathfrom IPython.display import display# Cargar datosdatos = pd.read_csv('datos/datos_sol_sus_acolite.csv')# Detectar columnas de bandasbandas = [col for col in datos.columns if col.startswith("B")]# Lista para guardar resultadosresultados = []# Parámetros para organización de gráficosn_bandas =len(bandas)ncols =3# Número de columnas de la grillanrows = math.ceil(n_bandas / ncols) # Calculamos cuántas filas se necesitan# Crear figurafig, axs = plt.subplots(nrows=nrows, ncols=ncols, figsize=(12, 4* nrows))axs = axs.flatten() # Asegura que podamos indexarlos como una listafor i, banda inenumerate(bandas):# Variables X = datos[[banda]] y = datos["sol_sus"]# División entrenamiento/prueba X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=False)# Ajuste del modelo modelo = LinearRegression().fit(X_train, y_train) y_train_pred = modelo.predict(X_train)# Métricas r2 = modelo.score(X_train, y_train) n =len(y_train) p = X_train.shape[1] r2_adj =1- (1- r2) * (n -1) / (n - p -1) rmse = np.sqrt(mean_squared_error(y_train, y_train_pred)) resultados.append({"Banda": banda,"R²": round(r2, 4),"R²_ajustado": round(r2_adj, 4),"RMSE": round(rmse, 4) })# Gráfico de entrenamiento ax = axs[i] ax.scatter(X_train, y_train, alpha=0.6, color="#9D50A6", label="Entrenamiento") ax.plot(X_train, y_train_pred, color="#17A77E", linewidth=1.8, label="Modelo") ax.set_title(f'{banda}\nR²={r2:.2f}', fontsize=10) ax.set_xlabel('Reflectancia', fontsize=8) ax.set_ylabel('Sol_Sus', fontsize=8) ax.tick_params(axis='both', which='major', labelsize=8) ax.legend(fontsize=7) ax.grid(True)# Eliminar ejes sobrantesfor j inrange(i +1, len(axs)): fig.delaxes(axs[j])plt.suptitle("Regresiones lineales por banda (entrenamiento)", fontsize=14)plt.tight_layout(rect=[0, 0, 1, 0.96])plt.show()# Tabla resumendf_resultados = ( pd.DataFrame(resultados) .round(3) .sort_values("R²", ascending=False) # la mejor banda queda arriba .reset_index(drop=True))def style_metrics(df: pd.DataFrame, caption: str): sty = ( df.style .format(precision=3) )# Ocultar índice según versión de pandasifhasattr(sty, "hide"): sty = sty.hide(axis="index")else: sty = sty.hide_index()# Agregar título y estilo del título sty = sty.set_caption(caption) sty = sty.set_table_styles([{"selector": "caption","props": [ ("caption-side", "top"), ("font-weight", "bold"), ("margin-bottom", "0.5em") ] }], overwrite=False)return sty# Mostrar tabla display(style_metrics(df_resultados,"Resumen de métricas por banda (entrenamiento)"))
Tabla 3: Resumen de métricas por banda (entrenamiento)
Banda
R²
R²_ajustado
RMSE
B05
0.183
0.167
33.457
B06
0.180
0.164
33.523
B07
0.173
0.157
33.651
B08
0.146
0.130
34.191
B8A
0.135
0.118
34.429
B04
0.094
0.076
35.230
B12
0.038
0.019
36.303
B11
0.033
0.015
36.389
B01
0.008
-0.011
36.865
B03
0.005
-0.014
36.917
B02
0.001
-0.018
36.985
Conclusiones del análisis por banda
Los resultados muestran que ninguna banda espectral logra predecir con alta precisión los sólidos suspendidos mediante un modelo lineal simple.
Las bandas B05, B06 y B07 presentan los mejores valores de R² (alrededor de 0.18), lo que indica que explican menos del 20 % de la variabilidad en los datos. El resto de las bandas tienen un desempeño aún menor, con errores (RMSE) relativamente altos en comparación con los valores observados.
Esto sugiere que la relación entre reflectancia y sólidos suspendidos no es lineal en escala natural y que será necesario explorar modelos más complejos o transformaciones de datos para mejorar la predicción.
2.5 Análisis de regresión lineal con transformación logarítmica
En esta etapa aplico una transformación logarítmica natural a las variables de reflectancia y sólidos suspendidos antes de ajustar los modelos de regresión lineal. Esta transformación permite estabilizar la varianza, linealizar relaciones no lineales y reducir el impacto de valores extremos.
El procedimiento es similar al análisis anterior, pero antes de entrenar el modelo se aplica log(x) a las columnas correspondientes. Para evitar problemas con valores cero, estos se reemplazan por NaN y se excluyen del análisis. Posteriormente, entreno modelos de regresión lineal simple por banda usando las variables transformadas.
El desempeño de cada modelo se evalúa mediante métricas en escala logarítmica (R²_log, RMSE_log) y en escala original, facilitando la comparación y la interpretación de resultados.
Código
import pandas as pdimport numpy as npimport matplotlib.pyplot as pltfrom sklearn.model_selection import train_test_splitfrom sklearn.linear_model import LinearRegressionfrom sklearn.metrics import mean_squared_error, r2_scoreimport mathimport seaborn as snsfrom IPython.display import display# 1. ── Cargar datos ────────────────────────────────────────────────datos = pd.read_csv("datos/datos_sol_sus_acolite.csv")bandas = [c for c in datos.columns if c.startswith("B")]# 2. ── Filtrar positivos (requisito para log) ──────────────────────mask_pos = (datos["sol_sus"] >0) & (datos[bandas] >0).all(axis=1)datos_log = datos.loc[mask_pos].copy()# 3. ── Transformar a log-natural ──────────────────────────────────datos_log["log_sol_sus"] = np.log(datos_log["sol_sus"])for b in bandas: datos_log[f"log_{b}"] = np.log(datos_log[b])log_bandas = [f"log_{b}"for b in bandas]# 4. ── Train / test split (descartamos test) ───────────────────────X_total = datos_log[log_bandas]y_total = datos_log["log_sol_sus"]X_train, _, y_train, _ = train_test_split( X_total, y_total, test_size=0.2, shuffle=False)# 5. ── Funciones de métricas ───────────────────────────────────────def metrics_log(model, X, y): y_pred = model.predict(X) rmse = np.sqrt(mean_squared_error(y, y_pred)) r2 = r2_score(y, y_pred) n, p = X.shape r2_adj =1- ((1- r2) * (n -1)) / (n - p -1)return rmse, r2, r2_adjdef metrics_original(model, X, y): y_pred_log = model.predict(X).ravel() y_pred_orig = np.exp(y_pred_log) y_true_orig = np.exp(y.values) rmse = np.sqrt(mean_squared_error(y_true_orig, y_pred_orig)) r2 = r2_score(y_true_orig, y_pred_orig) n, p = X.shape r2_adj =1- ((1- r2) * (n -1)) / (n - p -1)return rmse, r2, r2_adj# 6. ── Entrenar un modelo por banda ────────────────────────────────resultados, resultados_nolog = [], []lineas_log, lineas_nolog = {}, {}for b, lb inzip(bandas, log_bandas): Xtr = X_train[[lb]] reg = LinearRegression().fit(Xtr, y_train)# Métricas en log rmse_l, r2_l, r2a_l = metrics_log(reg, Xtr, y_train) resultados.append({"Banda": b, "R²_log": r2_l, "R²aj_log": r2a_l, "RMSE_log": rmse_l})# Métricas en original rmse_o, r2_o, r2a_o = metrics_original(reg, Xtr, y_train) resultados_nolog.append({"Banda": b,"R²": r2_o, "R²aj": r2a_o, "RMSE": rmse_o})# Rectas ordenadas order = np.argsort(Xtr.values.ravel()) x_log_ord = Xtr.values.ravel()[order] y_pred_log = reg.predict(Xtr).ravel()[order] lineas_log[b] = (x_log_ord, y_pred_log) x_orig_ord = np.exp(x_log_ord) y_pred_orig = np.exp(y_pred_log) lineas_nolog[b] = (x_orig_ord, y_pred_orig)# 7. ── DataFrames de métricas ──────────────────────────────────────res_df_log = (pd.DataFrame(resultados) .round(3) .sort_values("RMSE_log") .reset_index(drop=True))res_df_nolog = (pd.DataFrame(resultados_nolog) .round(3) .sort_values("RMSE") .reset_index(drop=True))# 8. ── Gráficos en escala original ────────────────────────────────n_cols, n_rows =3, math.ceil(len(bandas) /3)fig, axes = plt.subplots(n_rows, n_cols, figsize=(15, 4* n_rows))axes = axes.flatten()for idx, b inenumerate(bandas): ax = axes[idx] x_orig = np.exp(X_train[f"log_{b}"]) y_orig = np.exp(y_train) ax.scatter(x_orig, y_orig, alpha=0.6, color="#9D50A6", label="Entrenamiento") ax.plot(*lineas_nolog[b], lw=3, color="#17A77E", label="Modelo") ax.set(xlabel=b, ylabel="sol_sus", title=f"Banda {b} (escala original)") ax.legend() ax.grid(True)for j inrange(len(bandas), len(axes)): fig.delaxes(axes[j])plt.tight_layout()plt.show()# 9. ── Función para estilo de tablas ───────────────────────────────def style_metrics_sin_gradiente(df: pd.DataFrame, caption: str): sty = (df.style .format(precision=3)# .background_gradient(...) ← quitamos el gradiente .hide(axis="index") .set_caption(caption) .set_table_styles([{"selector": "caption","props": [("caption-side", "top"), ("font-weight", "bold"), ("margin-bottom","0.5em")] }], overwrite=False) )return sty# 10 ── Mostrar tablas con el nuevo estilo pastel ───────────────────────────────display(style_metrics_sin_gradiente(res_df_log,"Métricas en escala log‑log (ordenadas por RMSE_log)"))display(style_metrics_sin_gradiente(res_df_nolog,"Métricas en escala original (ordenadas por RMSE)"))
Tabla 4: Métricas en escala log‑log (ordenadas por RMSE_log)
Banda
R²_log
R²aj_log
RMSE_log
B06
0.304
0.281
0.318
B07
0.303
0.280
0.319
B05
0.272
0.248
0.326
B08
0.267
0.243
0.327
B8A
0.262
0.237
0.328
B04
0.166
0.138
0.349
B12
0.111
0.082
0.360
B11
0.065
0.034
0.369
B03
0.017
-0.016
0.378
B01
0.015
-0.018
0.379
B02
0.007
-0.026
0.380
Tabla 5: Métricas en escala original (ordenadas por RMSE)
Banda
R²
R²aj
RMSE
B06
0.198
0.171
31.661
B07
0.186
0.159
31.897
B05
0.178
0.150
32.049
B8A
0.156
0.128
32.472
B08
0.151
0.123
32.567
B04
0.077
0.046
33.962
B12
0.076
0.045
33.986
B11
0.043
0.011
34.576
B01
-0.028
-0.063
35.843
B03
-0.034
-0.068
35.940
B02
-0.035
-0.070
35.969
📄 Nota técnica
La transformación log(x) se usa frecuentemente cuando la relación entre variables es multiplicativa o cuando hay asimetría en la distribución.
Es fundamental reemplazar ceros por NaN antes de aplicar el logaritmo, ya que log(0) no está definido.
Las métricas obtenidas (R² y RMSE) se interpretan en la escala logarítmica y no son directamente comparables con las obtenidas en la escala original.
El R² ajustado penaliza la complejidad del modelo y ayuda a evaluar si el ajuste mejora más allá de lo que se esperaría por azar.
Conclusión del ajuste log–log por banda
Los resultados muestran que al aplicar la transformación logarítmica, las bandas B06, B07 y B05 obtienen los mejores ajustes, con valores de R²_log cercanos a 0.30 y errores (RMSE_log) más bajos comparados con las demás bandas. Esto indica una mejora en la capacidad predictiva respecto al análisis en escala natural.
Sin embargo, al evaluar las métricas en escala original, el mejor desempeño alcanza un R² máximo cercano a 0.20, lo que aún refleja una capacidad limitada para explicar la variabilidad de los sólidos suspendidos utilizando bandas individuales.
Estos resultados sugieren que la transformación logarítmica mejora el modelo, pero que para obtener predicciones más precisas será necesario explorar combinaciones de bandas, índices espectrales o modelos más complejos.
2.6 Selección de bandas mediante regresión logarítmica y AIC
Para determinar cuáles bandas espectrales aportan más información para predecir sólidos suspendidos, ajusté modelos de regresión lineal simple en escala logarítmica para cada banda individual. Además de las métricas habituales como el coeficiente de determinación (R²) y el error cuadrático medio (RMSE), utilicé el Criterio de Información de Akaike (AIC).
El AIC es una medida que evalúa la calidad del modelo penalizando la complejidad: un valor menor indica un mejor equilibrio entre ajuste y simplicidad. Esto ayuda a evitar sobreajuste al no favorecer modelos que mejoran el ajuste sólo por aumentar el número de parámetros sin aportar realmente información útil. Así, el AIC es una herramienta clave para seleccionar las bandas que efectivamente mejoran la predicción sin agregar complejidad innecesaria.
Código
import pandas as pdimport numpy as npfrom sklearn.model_selection import train_test_splitfrom sklearn.linear_model import LinearRegressionfrom sklearn.metrics import mean_squared_error, r2_scoreimport seaborn as snsfrom IPython.display import display# 1. Cargar datos y detectar bandasdatos = pd.read_csv("datos/datos_sol_sus_acolite.csv")bandas = [c for c in datos.columns if c.startswith("B")]# 2. Filtrar casos con valores positivos (requisito de log)mask_pos = (datos["sol_sus"] >0) & (datos[bandas] >0).all(axis=1)datos_log = datos.loc[mask_pos].copy()# 3. Transformar a log-naturaldatos_log["log_sol_sus"] = np.log(datos_log["sol_sus"])for b in bandas: datos_log[f"log_{b}"] = np.log(datos_log[b])log_bandas = [f"log_{b}"for b in bandas]# 4. Partición entrenamiento / test X_total = datos_log[log_bandas]y_total = datos_log["log_sol_sus"]X_train, _, y_train, _ = train_test_split( X_total, y_total, test_size=0.20, shuffle=False)# 5. Funciones auxiliares (AIC + métricas)def aic_gauss(model, X, y):""" AIC para regresión lineal con ruido ~ N(0, σ²) AIC = n · ln(RSS/n) + 2k """ n =len(y) k = X.shape[1] +1# predictores + intercepto rss = np.sum((y - model.predict(X)) **2)return n * np.log(rss / n) +2* kdef metrics_log(model, X, y): y_pred = model.predict(X) rmse = np.sqrt(mean_squared_error(y, y_pred)) r2 = r2_score(y, y_pred) n, p = X.shape r2_adj =1- ((1- r2) * (n -1)) / (n - p -1)return rmse, r2, r2_adj# 6. Ajuste univariable banda-a-banda + AICresultados = []for b in bandas: lb =f"log_{b}" Xtr = X_train[[lb]] reg = LinearRegression().fit(Xtr, y_train) rmse_log, r2_log, r2a_log = metrics_log(reg, Xtr, y_train) aic_val = aic_gauss(reg, Xtr, y_train) resultados.append({"Banda" : b,"R²" : r2_log,"R²aj" : r2a_log,"RMSE" : rmse_log,"AIC" : aic_val })# 7. Tabla ordenada por AICdf_log = (pd.DataFrame(resultados) .round(3) .sort_values("AIC") .reset_index(drop=True))def style_metrics_sin_gradiente(df: pd.DataFrame, caption: str): sty = (df.style .format(precision=3) .hide(axis="index") # oculta índice .set_caption(caption) .set_table_styles([{"selector": "caption","props": [("caption-side", "top"), ("font-weight", "bold"), ("margin-bottom", "0.5em")] }], overwrite=False) )return sty# 9 ── Mostrar tabla con el nuevo estilo pastel AIC ───────────────────────────display(style_metrics_sin_gradiente( df_log, "AIC y métricas en escala log‑log (ordenadas por AIC)"))
Tabla 6: AIC y métricas en escala log‑log (ordenadas por AIC)
Banda
R²
R²aj
RMSE
AIC
B06
0.304
0.281
0.318
-69.242
B07
0.303
0.280
0.319
-69.191
B05
0.272
0.248
0.326
-67.798
B08
0.267
0.243
0.327
-67.584
B8A
0.262
0.237
0.328
-67.355
B04
0.166
0.138
0.349
-63.448
B12
0.111
0.082
0.360
-61.415
B11
0.065
0.034
0.369
-59.781
B03
0.017
-0.016
0.378
-58.184
B01
0.015
-0.018
0.379
-58.112
B02
0.007
-0.026
0.380
-57.878
Conclusión sobre la selección de bandas con AIC
Los valores obtenidos muestran que las bandas B06 y B07 son las más relevantes para predecir los sólidos suspendidos, ya que presentan los valores más bajos de AIC (-69.24 y -69.19 respectivamente), indicando el mejor equilibrio entre calidad de ajuste y simplicidad del modelo.
Estas bandas también exhiben los mejores coeficientes de determinación (R²_log alrededor de 0.30) y los menores errores (RMSE_log cerca de 0.32).
El resto de las bandas tienen valores de AIC menos favorables, lo que sugiere que, individualmente, aportan menos información útil y no justifican la complejidad del modelo. Esto refuerza la idea de que para mejorar las predicciones será necesario combinar varias bandas o explorar modelos más complejos.
2.6.1 Selección de banda espectral con validación cruzada
En esta sección selecciono la banda espectral más adecuada para estimar los sólidos suspendidos (sol_sus) usando un modelo de regresión lineal simple en escala logarítmica. A diferencia de análisis anteriores, ahora incorporo una validación cruzada que respeta la estructura temporal de los datos, para asegurar que el modelo sea confiable al aplicarse a nuevas campañas de muestreo.
¿Por qué validación Leave-One-Group-Out (LOGO)?
Los datos provienen de campañas realizadas en diferentes fechas (fecha). Esto genera grupos independientes, ya que las condiciones del río pueden variar entre campañas. Usar una validación cruzada tradicional como K-Fold o Shuffle Split mezclaría muestras de una misma campaña entre entrenamiento y prueba, lo que daría una estimación demasiado optimista del desempeño.
En cambio, LeaveOneGroupOut (LOGO) deja una fecha completa afuera en cada iteración, lo que permite evaluar si el modelo puede generalizar a campañas no vistas. Esta estrategia evita filtraciones de información (data leakage) entre los datos de entrenamiento y validación, y resulta mucho más realista en contextos ambientales, donde el objetivo es predecir nuevas campañas a partir de mediciones anteriores.
Descripción del procedimiento:
Preprocesamiento: aplico log() a las bandas y a la variable objetivo (sol_sus) y elimino ceros, para linearizar y estabilizar la relación.
Modelado: ajusto un modelo de regresión lineal simple por cada banda, normalizando los datos con StandardScaler (pipeline). Para cada modelo calculo las métricas de validación cruzada: RMSE, R², R² ajustado y el AIC, que penaliza la complejidad del modelo.
Selección: ordeno las bandas por AIC y RMSE, y elijo la mejor.
Evaluación final:
Estimo el desempeño real sobre la campaña más reciente (test).
Realizo un análisis bootstrap (500 réplicas) para obtener intervalos de confianza de los parámetros y métricas.
Muestro gráficos de predicho vs observado en escala log y original.
Este enfoque no solo me permite elegir la banda más informativa, sino que también evalúa la estabilidad del modelo y su capacidad de generalizar en el tiempo, algo clave para el monitoreo operativo de calidad de agua.
Selecciono la banda B07 como la mejor candidata para estimar los sólidos suspendidos. Esta elección se basa en el valor más bajo del AIC y el menor error promedio en validación cruzada (RMSE_CV). El modelo resultante indica que, a mayor reflectancia en la banda B07, tiende a haber mayor concentración estimada de sólidos suspendidos.
Sin embargo, al evaluar el modelo con los datos más recientes (campaña de prueba), el desempeño cae considerablemente. El R² es negativo, lo que significa que el modelo predice peor que una estimación promedio constante. En otras palabras, el modelo no generaliza bien a nuevos datos.
Para tener una idea más completa de su estabilidad, aplico la técnica de bootstrap con 500 repeticiones. Esto muestra que si bien el coeficiente de la banda B07 es consistentemente positivo (lo que confirma una relación), el rendimiento general del modelo es bajo y muy variable.
Concluyo que, aunque B07 es la mejor banda individual, no alcanza para predecir sólidos suspendidos con precisión. Para mejorar el modelo, es necesario combinar varias bandas, usar índices espectrales o aplicar modelos más complejos.
2.7 Selección de bandas con Forward Selection basado en AIC
En esta sección realizo la selección de un conjunto de bandas espectrales para predecir la concentración de sólidos suspendidos (sol_sus) mediante un modelo de regresión lineal múltiple en escala logarítmica. Utilizo el método de Forward Selection guiado por el criterio de información de Akaike (AIC), que permite incorporar progresivamente aquellas bandas cuya inclusión mejora significativamente el ajuste del modelo sin incrementar de forma innecesaria su complejidad.
En esta sección selecciono un conjunto de bandas espectrales para predecir los sólidos suspendidos (sol_sus) mediante regresión lineal múltiple en escala logarítmica. Utilizo el método de Forward Selection, que comienza con un modelo vacío e incorpora progresivamente las bandas que más reducen el AIC, siempre que la mejora sea significativa (ΔAIC ≥ 1).
Transformo todas las variables a escala logarítmica para linealizar relaciones y reducir la variabilidad. Reservo la última fecha como conjunto de prueba, y entreno el modelo con el resto de los datos. Durante la selección, evalúo el desempeño con validación cruzada agrupada (GroupKFold), que respeta las campañas de muestreo y permite estimar la capacidad del modelo para generalizar a nuevas fechas.
Finalmente, ajusto el modelo con las bandas seleccionadas, calculo métricas en train y test, y realizo un análisis bootstrap para obtener intervalos de confianza de los coeficientes. Esto me permite construir un modelo interpretable, robusto y ajustado a las condiciones reales del monitoreo ambiental.
Selecciono las bandas B07 y B02 como las mejores para estimar los sólidos suspendidos, basándome en la reducción significativa del AIC y la mejora en el error promedio de validación cruzada (RMSE_CV). El modelo indica que la reflectancia en B07 tiene un efecto positivo en la concentración estimada, mientras que B02 aporta un efecto negativo, lo que puede relacionarse con propiedades ópticas del agua y sedimentos.
Al evaluar el modelo con datos nuevos (campaña de prueba), el desempeño es razonablemente bueno, con un R² positivo que refleja una capacidad aceptable de generalización. El análisis bootstrap confirma la estabilidad de los coeficientes y las métricas, entregando confianza en la robustez del modelo.
En resumen, combinar estas dos bandas mejora sustancialmente la predicción en comparación con usar una sola banda, pero para lograr mayor precisión probablemente sea necesario explorar índices espectrales o modelos más avanzados.
2.8 Selección de variables mediante Forward Selection con AIC usando razones simples de bandas
En esta etapa exploro combinaciones de variables formadas por razones simples entre bandas espectrales, todas en escala logarítmica, para mejorar la predicción de los sólidos suspendidos (sol_sus). Aplico un método de Forward Selection guiado por el criterio de información de Akaike (AIC), que me permite construir un modelo parsimonioso incorporando solo aquellas razones que aportan mejoras significativas en el ajuste.
El análisis se realiza sobre el conjunto de entrenamiento definido previamente, manteniendo la partición temporal para test y respetando la estructura grupal en la validación cruzada. Además, implemento un análisis bootstrap para estimar la incertidumbre en los coeficientes y evaluar la estabilidad del modelo.
Con este procedimiento busco identificar relaciones espectrales más complejas que mejor expliquen la variabilidad de sólidos suspendidos, maximizando la capacidad predictiva sin sobreajustar, para un monitoreo ambiental más preciso y confiable.
En esta etapa, al usar razones simples de bandas espectrales en lugar de bandas individuales, logro un modelo con mejor ajuste en entrenamiento, con un R² de 0.784 frente a 0.766 del modelo previo. También el error promedio en la validación cruzada mejora ligeramente (RMSE CV 0.119 vs. 0.125), lo que indica una selección más eficiente de variables.
Sin embargo, al evaluar con los datos de la campaña más reciente, el desempeño baja significativamente: el R² cae a 0.212, mucho menor que el 0.484 obtenido anteriormente, y el RMSE aumenta, mostrando que la capacidad de generalización sigue siendo limitada.
El análisis bootstrap confirma que la razón log(B02/B07) tiene un coeficiente negativo estable y significativo, mientras que log(B05/B11) aporta muy poco y su coeficiente incluye cero en el intervalo de confianza, indicando que su efecto es marginal.
Desde el punto de vista físico, el modelo indica que a medida que aumenta la razón B02/B07, la concentración de sólidos suspendidos tiende a disminuir, mientras que la razón B05/B11 tiene un efecto muy débil e incierto.
A pesar de que el uso de razones simples resulta en un modelo más compacto y con una mejora modesta en el ajuste, la capacidad para predecir sólidos suspendidos en campañas futuras continúa siendo limitada.
2.8.1 Ajuste del umbral ΔAIC para selección de variables
Para evaluar la robustez y simplicidad del modelo, aumento el umbral mínimo de mejora en el criterio AIC (ΔAIC) de 2 a 4. Esto implica que la selección hacia adelante solo incorpora nuevas variables si la reducción en AIC es al menos 4, promoviendo modelos más parsimoniosos y evitando incluir variables con mejoras marginales que podrían no generalizar bien.
Este ajuste busca balancear mejor el compromiso entre ajuste y complejidad, reduciendo el riesgo de sobreajuste, especialmente relevante al aplicar el modelo a campañas de muestreo futuras.
TRAIN campañas : 5
TEST campaña : 2025-06-09 (12 muestras)
Tabla 16: Historial de selección
Paso
Bandas
AIC
RMSE
R²
R²aj
1
log(B02/B07)
-106.889
0.148
0.739
0.729
Bandas seleccionadas:log(B02/B07)
GroupKFold CV (log) RMSE: 0.141 ± 0.053
Tabla 17: Métricas finales
Métrica
Train log
Test log
Train orig
Test orig
RMSE
0.148
0.391
9.707
46.109
R²
0.739
0.202
0.797
0.013
Tabla 18: Bootstrap (coeficientes, 95 % CI)
Banda
Media β
IC 95% inf.
IC 95% sup.
log(B02/B07)
-0.5622
-0.6158
-0.5001
Ecuación (log)
log(sol_sus) = 4.2877 + -0.5625·log(B02/B07)
Ecuación (original)
sol_sus = 72.8024 · log(B02/B07)^-0.5625
Impacto del aumento del umbral ΔAIC en la selección de variables
Al aumentar ΔAIC a 4, el modelo se simplifica seleccionando solo la razón log(B02/B07), manteniendo un buen ajuste en entrenamiento (R² = 0.739). Sin embargo, el desempeño en la campaña de prueba es bajo, con R² cercano a cero, indicando poca capacidad de generalización. El análisis bootstrap muestra estabilidad en el coeficiente seleccionado.
Por lo tanto, un mayor ΔAIC favorece modelos más simples, pero no mejora la predicción en nuevas campañas, por lo que es necesario explorar modelos más complejos o incluir más variables.
2.9 Selección de bandas, razones y potencias para modelar sólidos suspendidos
En esta sección amplío el conjunto de variables predictoras al incluir no solo logaritmos simples de bandas espectrales, sino también sus potencias cuadráticas y cúbicas, así como razones logarítmicas entre bandas. El objetivo es capturar relaciones no lineales y combinadas que podrían mejorar la capacidad predictiva del modelo.
Utilizo una estrategia de selección hacia adelante basada en el criterio AIC (ΔAIC ≥ 4), priorizando modelos más parsimoniosos pero informativos. La validación se realiza con partición temporal (hold-out) y GroupKFold, y se evalúa la estabilidad del modelo mediante bootstrap de coeficientes (IC del 95 %). Finalmente, se comparan las predicciones con observaciones en escalas logarítmica y original.
Índice de agua → min: -12.426, max: 17.432
% de píxeles clasificados como agua: 76.45%
--- Umbral: 0 ---
Estadísticas del índice de agua:
Min: -12.4261, Max: 17.4322, Media: 0.6814
Estadísticas de la máscara de agua:
Valores únicos: [0 1] (1=agua, 0=tierra, 255=nodata)
Estadísticas del parámetro sólidos suspendidos:
Min: 19.9269, Max: 51112.1016, Media: 139.5654
Índice de agua → min: -12.426, max: 17.432
% de píxeles clasificados como agua: 74.87%
Índice de agua → min: -12.426, max: 17.432
% de píxeles clasificados como agua: 74.74%
Índice de agua → min: -12.426, max: 17.432
% de píxeles clasificados como agua: 74.73%
Índice de agua → min: -12.426, max: 17.432
% de píxeles clasificados como agua: 74.71%
Referencias
[1]
D. C. R. Ramírez, «Método de Estimación de Sólidos Suspendidos Totales como Indicador de la Calidad del Agua Mediante Imágenes Satelitales». Universidad Nacional de Colombia, Facultad de Ciencias Agrarias, 2017.
[2]
G. D. J. L., J. Stephan, y Delance-Martinic., «Determinación del parámetro sólidos suspendidos totales (sst) mediante imágenes de sensores ópticos en un tramo de la cuenca media del río Bogotá (Colombia).», UD y la Geomática, pp. 19-27, 2014, doi: 10.14483/udistrital.jour.udgeo.2014.9.a02.
[3]
A. Cruz-Retana, C. Fonseca-Ortiz, R. Becerril-Piña, M. A. Gómez-Albores, M. Hernández-Téllez, y R. Arévalo-Mejia, «Characterization of spectral reflectance and TSS concentration in continental water bodies worldwide», vol. 1. pp. 4-18, 2023.
[4]
E. E. C. Vargas, «Modelamiento de Relaciones entre Parámetros Fisicoquímicos y Microbiológicos en Aguas de la Bahía Interior del Lago Titicaca-Puno (Perú) mediante Árboles de Predicción», Revista Tecnica De La Facultad De Ingenieria Universidad Del Zulia, vol. 44, pp. 154-168, ago. 2021, doi: 10.22209/rt.v44n3a02.
[5]
M. Moeini, A. Shojaeizadeh, y M. Geza, «Supervised machine learning for estimation of total suspended solids in urban watersheds», Water (Switzerland), vol. 13, ene. 2021, doi: 10.3390/w13020147.
[6]
L. S. Kupssinskü, T. T. Guimarães, E. M. D. Souza, y D. C. Zanotta, «A method for chlorophyll-a and suspended solids prediction through remote sensing and machine learning», Sensors (Switzerland), vol. 20, abr. 2020, doi: 10.3390/s20072125.